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' We investigate the stationary states of one-dimensional driven diffusive systems, coupled to 

boundary reservoirs with fixed particle densities. We argue that the generic phase diagram is 
. governed by an extremal principle for the macroscopic current irrespective of the local dynamics. 

In particular, we predict a minimal current phase for systems with local minimum in the current- 
density relation. This phase is explained by a dynamical phenomenon, the branching and coalescence 
of shocks, Monte-Carlo simulations confirm the theoretical scenario. 



A recurrent problem in the investigation of many-body systems far from equilibrium is posed by the coupling of 
O ' a driven particle system with locally conserved particle number to external reservoirs with which the system can 
exchange particles at its boundaries. In the presence of a driving force a particle current will be maintained and 
hence the system will always remain in an nonequilibrium stationary state characterized by some bulk density and 
the corresponding particle current. While for periodic boundaries the density is a fixed quantity, the experimentally 
more relevant scenario of open boundaries naturally leads to the question of steady-state selection, i.e. the question 



, which stationary bulk density the system will assume as a function of the boundary densities [pj . In the topologically 



C3 



I 



simplest case of quasi one-dimensional systems this is of importance for the understanding of many-body systems in 
which the dynamic degrees of freedom reduce to effectively one dimension as e.g. in traffic flow ||, kinetics of protein 
synthesis || , or diffusion in narrow channels Q] . 

Within a phenomenological approach this problem was first addressed in general terms by Krug J5| who postulated 
£h ' a maximal-current principle for the specific case where the density p + at the right boundary to which particles 
are driven is kept at zero. The exact solution of the totally asymmetric simple exclusion process (TASEP) || for 
arbitrary left and right boundary densities p~ ,p + confirms and extends the results by Krug. The complete phase 
diagram comprises boundary-induced phase transitions of first order between a low-density phase and a high density 
phase, and a second-order transition from these phases to a maximal current phase [p|dlo|. Analysis of the density 
profile j^] provides insight into the dynamical mechanisms that lead to these phase transitions and shows that the phase 
diagram is generic for systems with a single maximum in the current-density relation [fllf . Experimental evidence for 
, the first-order transition is found in the process of biopolymerization for which the TASEP with open boundaries was 
04 ■ originally invented as a simple model |l2j , and more directly also in recent measurements of highway traffic close to 
on-ramps [ |p3|]b"f| . Renormalization group studies |l5|] indicate universality of the second order phase transition. 

In this Letter we develop a dynamical approach to generic driven one-component systems with several maxima in 
the current-density relation, and we show that a novel phase of rather unexpected nature appears: For a certain range 
of boundary densities the steady state carries the minimal current between two maxima even though both boundary 
densities support a higher current. We shall refer to this phase as minimal current phase. More generally, we claim 
that the current always obeys the extremal principle 

x3 ; 

j = max j(p) for p > p + (1) 
O . P e[p+, P -] 
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j= min j(p) for p < p + . (2) 

pe[p~,p+] 



To understand the origin of this extremal principle we first note that in the absence of detailed balance stationary 
behavior cannot be understood in terms of a free energy, but has to be derived from the system dynamics. For 
definiteness consider a driven lattice gas with hard-core repulsion. At p = 1 no hopping can take place and hence the 
current vanishes. Two local maxima can arise as the result of sufficiently strong repulsion between nearest neighbor 
particles as opposed to the pure on-site repulsion of the usual TASEP which leads to a single maximum (Fig. |l|). 

At first sight one might not expect such a little change in the interaction radius of the particles to affect the phase 
diagram. However, the extremal principle (]l|), (||) predicts that the full phase diagram (Fig. ||) generically consists 
of seven distinct phases, including two maximal current phases with bulk densities corresponding to the respective 
maxima of the current Jl6| and the minimal current phase in a regime defined by 

j(P + ),j(P~) > j(Prmn); P~ < Pmin < P + ■ (3) 
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Here the system organizes itself into a state with bulk density pbuik corresponding to the local minimum of the current. 
As in the maximal current phases no fine-tuning of the boundary densities is required. 

In order to understand the basic mechanisms which determine the steady state selection we first follow Ref. @ 
and consider the collective velocity 

v c =j'(p) (4) 

which is the velocity of the center of mass of a local perturbation in a homogeneous, stationary background of density 
p (Fig. 3(a)). A second quantity of interest is the shock velocity 

h ~ h /,x 
v s = (5) 

P2 - Pi 

of a 'domain wall' between two stationary regions of densities p±^ (Fig. 3(b)). Notice that both velocities may change 
sign, depending on the densities. 

These velocities are sufficient to understand the phase diagram of systems with a single maximum in the current 
(□J. Further to these notions we need to introduce here the concepts of coalescence and branching of shocks. A single 
large shock (with a large density difference P2 — pi) may be understood as being composed of subsequent smaller 
shocks with narrow plateaux at each level of density (Fig. |^). In the usual asymmetric exclusion process these shocks 
travel with different relative speeds such that they coalesce and form a 'bound state' equivalent to a single shock Jl7| . 
In the present situation, however, the minimum in the current-density relation leads to more complicated dynamics. 
A closer investigation of Eq. (||) shows that depending on p\^ a single shock may branch into two distinct shocks, 
moving away from each other, discussed below. 

With these observations the dynamical origin of the phase transition lines can be understood by considering the 
time evolution of judiciously chosen initial states. Because of ergodicity, the steady state does not depend on the 
initial conditions and a specific choice involves no loss of generality. We turn our attention to a line p + = c with 
Pmin < c < pi, in the phase diagram which crosses the minimal current phase. Along this line it is convenient to 
consider an initial configuration with a shock with densities p~ and p + on the left and on the right respectively, which 
is composed of many narrow subsequent shocks at various levels of intermediate densities (Fig. ^). 

(i) Let us start with the point with equal boundary densities in which case the system evolves into a steady state 
with the same bulk density pbuik = P~ = P + ■ 

(ii) Lower p~ slightly below p + with just a single shock separating both regions. According to (||) the shock travels 
with speed v s = (j + —j~)/{p + — p~) > to the right, making the bulk density equal to p~~ . At the same time, small 
disturbances will according to (||) also drift to the right, as v c = j'(ft~) > in this region, thus stabilizing the single 
shock. 

(iii) Now, lower p~ below p m i n - While the shock velocity v s is still positive, so that one expects the shock to 
move to the right, the collective velocity v c = j'(p~) < indicates that disturbances will spread to the left. This 
discrepancy marks the failure of a single shock scenario. In order to resolve it, consider for simplicity p~ = pi and 
return to the picture with many subsequent shocks at each density level between p~ and p + (Figfl). Eq. (||) shows 
that small shocks below p m i n will move to the left, while those above p m in will move to the rightjH The leftmost of 
the left-moving shocks will merge in a single one, and so will the rightmost of the right-moving shocks. The result 
are two single shocks and moving different directions and thus expanding the region with the density pbuik = Pmin- 
The system is in the minimal current phase. This picture is well supported by the Monte-Carlo simulations shown 
in Fig. |5|, demonstrating the branching of a single shock into two distinct shocks moving in opposite directions. The 
minimal current phase will persist for any left boundary density in the range p~ G [pi , p m in] ■ Notice that the change 
of bulk density is continuous across the point p~ = p m i n to the minimal current phase, so the transition is of the 
second order. 

(iv) As we lower p~ below p\, the shock velocity v s — (j m in — j(p~))/(pmin ~ p) > becomes positive. The shock 
is moving to the right, leading to a low density phase with bulk density pbuik = p~ which drops discontinuously from 
Pbuik — Pmin at p~ = p% + to pbuik = P~ at p~ = p\ — 0. The system undergoes a first order phase transition. On 



1 Actually, the small shock between the levels p^,p\ will be stable only if v c (p^) < v s < Vc(p%), see Eq.([||j). Some of small 
shocks satisfy the above condition and some do not. This brings about additional structure to the resulting shock, which will 
be discussed elsewhere 
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the transition line the shock performs an unbiased random walk, separating coexisting regions of densities p m in and 
p~ respectively. 

(v) Let us start again from p~ = p + and now increase p~ . Until one reaches p~ = pi,-, the collective velocity 
v c = j'(p~) > is positive, leading to p bu i k = p~ ■ 

(vi) As soon as p~ crosses the point p~ = p%, the sign of collective velocity v c changes and the overfeeding effect 
H occurs: a perturbation from the left does not spread into the bulk and therefore further increase of the left 
boundary density does not increase the bulk density. The system enters the maximal current phase II through a 
second-order transition. 

Using analogous arguments one constructs the complete phase diagram (Fig. |^) and obtains the extremal principle 
fljj |). The velocities (g), (§) which deter mine the phase transition lines follow from the current-density relation. This 
behavior can be checked with Monte Carlo simulations. A model with two maxima of the current is a TASEP with 
next nearest neighbor interaction defined by the bulk hopping rates (see A26,A27 in [[l8]) 

0100^0010 with rate 1 + 5 (6) 

1100- ^1010 with rate 1 + e (7) 
0101->0011 with rate 1 - e (8) 

1101— >1011 with rate 1-5 (9) 

with |e| < 1; \S\ < 1. For sufficiently strong repulsive interaction 1 — e <C 1 the current at half- filling is strongly 
suppressed, which brings about two maxima structure in the current-density relation (Fig. |l|). The limit e = 1 leads to 
jmin = 0. For negative e (attractive interaction) the current-density relation always has a single maximum. The other 
parameter, 5, is responsible for the particle-hole symmetry. 5 — corresponds to symmetrical graph j(p) = j(l — p). 
(5^0 breaks the particle-hole symmetry in favor of larger particles current ( <5 > 0) or larger vacancies current ( 
5 < 0). In particular, 5 ^ is responsible for different heights of two maxima on Fig. nj. 

The injection at the left boundary site 1 and extraction of particles at the right boundary site L is chosen to 
correspond to coupling to boundary reservoirs with densities p respectively. Along the line p + = p~ the stationary 
distribution is then exactly given by the equilibrium distribution of an 1-D Ising model with boundary fields and the 
bulk field such that the density profile is constant with density p = p + = p~ (T^|. The current j = (1 + 5) (0100) + 
(1 + e)( 1100) + (1 — e)( 0101 ) + (1 — 5)( 1101 ) as a function of the density can be calculated exactly using standard 
transfer matrix techniques. The exact graph is shown in Fig. |l| for specific values of the hopping rates. 

We performed Monte Carlo simulations for systems of sizes L from 100 to 1000. Densities and currents were 
averaged over at least 50L rounds, and averaged over 7 different histories. We consider the bulk density pbuik as 
order parameter. In analogy to equilibrium phase transitions a singularity (jump discontinuity) in the first derivatives 
of Pbuik signals a first (second) order phase transition. In the finite system we determined the location of the first 
(second) order transition by the appearance of peak (jump) in the first derivatives of the bulk density pbuik(p + , p ) 
with respect to p + and p~ . As initial state we chose either the empty or the completely filled lattice, whichever gave 
the faster convervence. Despite finite size effects, precise analysis of which requires further investigations, the overall 
agreement of the simulated phase diagram with the predicted one is very good already for L = 150 (Fig.0)H. 

We conclude that the interplay of density fluctuations and shock diffusion, coalescence and branching resp. as 
described above represents the basic mechanisms which determine the steady state selection (0,0) of driven diffusive 
systems with a nonlinear current density relation. A surprising phenomenon is the occurrence of the self-organized 
minimal current phase. Since little reference is made to the precise nature of the dynamics we argue that the 
phase diagram is generic and hence knowledge of the macroscopic current-density relation of a given physical system 
is sufficient to calculate the exact nonequilibrium phase transition lines which determine the density of the bulk 
stationary state. For systems with more than two maxima in the current the interplay of more than two shocks has 
to be considered. 

A fundamentally different scenario seems likely in systems with long-range bulk correlations, which, however, is 
untypical in one dimension. Perhaps the most interesting open question which one might address in a similar fashion is 
the question of steady-state selection in systems with two or more local conservation laws. This could provide insight 
into the mechanisms responsible for boundary-induced spontaneous symmetry breaking |pof in 1-D driven diffusive 
systems with open boundaries. 



2 We have chosen the set of parameters e and 5, which shows distinct minimum between two maxima. For the other choices 
of e, 5, the phase diagram can be obtained in the same way, e.g. by current-density relation and extremal principles Eqs.dlpl). 
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FIG. 1. Exact current-density relation of the TASEP with nearest-neighbor interaction for e = 0.995, S — 0.2 (Eqs.M-(H)). 
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FIG. 2. Exact phase diagram as a function of the boundary densities p , (p + ). Full (bold) lines indicate phase transitions of 
second (first) order. Circles show the results of Monte-Carlo simulations of a system with 150 sites where e = 0.995, 5 = 0.2. 
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FIG. 3. (a) Diffusive spreading of a density perturbation in the stationary state at two times ti > ti. The collective velocity 
describes the motion of the center of mass of the perturbation, (b) Motion of a shock. To the left of the domain wall particles 
are distributed homogeneously with an average density pi . To the right of the domain wall the background density is P2 
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FIG. 4. Schematical drawing of the decomposition of a large shock into small shocks and their velocities, leading to branching 
and coalescence. Here p\ < p~ < p m in < P + < P*i- 
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FIG. 5. Snapshots of a particle density distribution at the initial moment of time and after 300 Monte-Carlo steps, showing 
expansion of the minimal current phase. Simulated is the system of 150 sites, with particles initially distributed with average 
density p~ — 0.1 (p + = 0.85) on the left (on the right). 3000 different histories are averaged over. 



6 



